Bifurcation analysis of the plane sheet pinch 
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A numerical bifurcation analysis of the electrically driven plane sheet pinch is presented. The 
electrical conductivity varies across the sheet such as to allow instability of the quiescent basic 
state at some critical Hartmann number. The most unstable perturbation is the two-dimensional 
tearing mode. Restricting the whole problem to two spatial dimensions, this mode is followed up to a 
time-asymptotic steady state, which proves to be sensitive to three-dimensional perturbations even 
close to the point where the primary instability sets in. A comprehensive three-dimensional stability 
analysis of the two-dimensional steady tearing-mode state is performed by varying parameters of 
the sheet pinch. The instability with respect to three-dimensional perturbations is suppressed by a 
sufficiently strong magnetic field in the invariant direction of the equilibrium. For a special choice 
of the system parameters, the unstably perturbed state is followed up in its nonlinear evolution and 
is found to approach a three-dimensional steady state. 



I. INTRODUCTION 

Finite-resistivity plasma instabilities play an impor- 
tant role for the release of stored magnetic energy in 
many astrophysical objects. They also restrict the plasma 
stability in several fusion devices The simplest con- 
figuration in which they appear is the plane sheet pinch. 
In a pinch a conducting fluid can be held together by 
the action of an electric current passing through it with 
the pressure gradients being balanced by the Lorentz 
force. Of special interest is the resistive tearing insta- 
bility, which was studied by Furth et al. Q by using 
a boundary layer approach and afterwards numerically 
without making the boundary-layer approximation ||-|| . 
All these studies refer to the infinite Hartmann number 
case because of their neglect of the kinematic viscosity. 
The Hartmann number Ha, which is the geometric mean 
of two Reynolds- like numbers, one being kinetic and the 
other magnetic, is the essential parameter that deter- 
mines the global stability boundaries of the plane sheet 
pinch as well as those of its cylindrical counterpart 
Thus kinematic viscosity has to be included. 

A recent sheet pinch study || has been done with spa- 
tially and temporally uniform kinematic viscosity and 
magnetic diffusivity, and with impenetrable stress-free 
boundaries. It is found that the quiescent ground state 
(in which the current density is uniform and the magnetic 
field profile across the sheet is linear) remains stable, no 
matter how strong the driving electric field. This study 
was extended to the case of magnetic diffusivity vary- 
ing across the sheet, which results in the profiles of the 
equilibrium magnetic field deviating from linear behav- 
ior. In particular, the conductivity profile may be chosen 
such that the magnetic- field and/or the current profile 
have inflection points. A Squire theorem could be proven 
for this configuration || whose stability depends on the 
Hartmann number, the degree of current concentration 
about the midplane of the sheet, and on the magnetic 
shear (i.e, the asymmetry of the equilibrium magnetic 



field) 0. 

A stability analysis can be considered as part of a bi- 
furcation analysis, which will be provided for the cases of 
two as well as of three spatial dimensions in the present 
paper. In a bifurcation analysis one tries to determine 
the set of possible time-asymptotic states, the attractors, 
for given values of the system parameters. The bifurca- 
tions from a static sheet-pinch equilibrium have previ- 
ously been studied for the case of two spatial dimensions 
||Tiyi2|| , Grauer (l3) investigated the interaction of two 
different tearing modes in the two-dimensional slab ge- 
ometry by reducing the dynamics at the bifurcation point 
to that on a four-dimensional center manifold. The new 
time-asymptotic states were found to be of the tearing- 
mode type, but e.g. also traveling waves were found. 

We note here that even though with increasing Hart- 
mann number the equilibrium becomes first unstable to 
two-dimensional perturbations (according to the Squire 
theorem), the new final states may be three-dimensional. 
This is one of the problems the present paper is ad- 
dressed to. The first question arising is which type of two- 
dimensional (2D) time-asymptotic state develops nonlin- 
early from the tearing mode when the whole problem 
is restricted to two spatial dimensions. The bifurcation 
studies p-|l3|] predict steady states for the generic cases. 
Similarly, large scale perturbations of a sheared magnetic 
field equilibrium were found to result in a final tearing- 
mode type stationary state via multiple coalescence of 
the magnetic island structures A question coming 
up then is whether the 2D time-asymptotic states are 
stable with respect to three-dimensional perturbations. 
If not, how do the stability properties depend on pa- 
rameters like the strength of a constant external mag- 
netic field or the wavelength of perturbations in the in- 
variant direction of the 2D state? Finally, what are the 
characteristic properties of the three-dimensional time- 
asymptotic states, when they manifest? In the present 
paper for the first time a comprehensive stability anal- 
ysis of the two-dimensional time-asymptotic states that 
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develop from the tearing mode is presented. For the case 
of a spatially uniform resistivity these problems were 
addressed in numerical studies of the magnetohydrody- 
namic (MHD) equations by Dahlburg et al. who 
found two-dimensional quasi-equilibria of the tearing- 
mode type to be unstable to three-dimensional pertur- 
bations. Secondary three-dimensional instabilities were 
similarly observed for non-static primary states with, in 
addition to a sheared magnetic field, a pressure-driven 
jet-like flow fl7| , and their nonlinear development was 
proposed as a scenario for the transition to MHD turbu- 
lence. 

The outline of the paper is as follows. In Sec. II the 
physical model is introduced. The MHD equations as well 
as boundary and initial conditions are discussed. In Sec. 
Ill we provide the results of the bifurcation analysis. First 
the 2D results are discussed. After that we investigate 
the linear stability of 2D time-asymptotic states with re- 
spect to 3D perturbations. For appropriate choices of the 
external parameters the 2D states prove to be 3D unsta- 
ble. The 3D asymptotics is investigated by a full three- 
dimensional long-time simulation of the pinch dynamics. 
Finally we discuss our results and end with an outlook 
in Sec. IV. 



and E are normalized by L, Bq, va, ta = L/va, Po^a, 
Bq/(hoL), and Bqva, respectively. The quantity E is the 
electric field. Equations (fil) and (0) then become 



— = -(v ■ V)v + M V v - Vp + J x B, 



f)Ti 

— = -Vx (S-^J-vxB), 



where 



M = ^andS=^ 
v ?7o 



(5) 



(6) 



(7) 



are Reynolds-like numbers based on the Alfven velocity: 
S is the Lundquist number and M its viscous analogue. 
The geometric mean of the two Reynolds-like numbers 
gives the Hartmann number, 



Ha = VMS. (8) 
Finally, the dimensionless Ohm's law becomes 

r 1 fJ=E + vxB. (9) 



II. PHYSICAL MODEL 



B. Boundary conditions and static equilibrium 



A. MHD equations 

We use the nonrelativistic, incompressible MHD equa- 
tions, 

p + (v • V)vj = pvV 2 v- Vp + J x B, (1) 



We use Cartesian coordinates X\, x%, x 3 and consider 
our magnetofluid in the slab < X\ < 1. Xi is referred 
to as the cross-sheet coordinate. In the x 2 and x% direc- 
tions periodic boundary conditions with periods L 2 and 
L3, respectively, are used. 

The boundary planes are assumed to be impenetrable 
and stress- free, i.e., 



<9B 

— = -V x (»7/i J - v x B), 



V-v = 0, V-B = 0, 



(2) 



(3) 



where v is the fluid velocity, B the magnetic induction, /z 
the magnetic permeability in a vacuum, J = Vx B//xq 
the electric current density, p the mass density, p the 
pressure, v the kinematic viscosity, and 77 the magnetic 
diffusivity. While p and v are assumed constant, r\ varies 
spatially: 



?7(x) = ?7o?Kx), 



(4) 



where 770 is a dimensional constant and 77 (x) a dimen- 
sionless function of position. 

Let the pinch width (sheet thickness) L = L\ and some 
yet arbitrary field strength Bq be used as the units of 
length and magnetic induction. Writing va — B / ^PoP 
for the Alfven velocity corresponding to Bq, we transform 
to dimensionless quantities. Specifically x, B, v, t, p, J, 



vi 



8V2 

dx\ 



dv 3 
dx\ 



= at x\ = 0, 1. 



(10) 



The system is driven by an electric field of strength E* 
in the X3 direction, which can be prescribed only on the 
boundary. We further assume that there is no magnetic 
flux through the boundary, 



B x =0 at a;i=0,l. 



(11) 



Conditions ( |l0| ) and (|Tl]) imply that the tangential 
components of v x B on the boundary planes vanish, 
so that according to Eq. @ 

E* S 

■h = 0, J 3 = at a:i = 0, 1, (12) 

where fji, is the value of 77 on the boundaries. The bound- 
ary conditions for the tangential components of B then 
become 



dB 2 E*S dB* 



dxi 1% ' dxi 



= at xi = 0, 1. 



(13) 
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A detailed discussion of these boundary conditions is 
found in Ref. §. 

Any stationary state with the fluid at rest has to satisfy 
the equations 



— Vp + J x B = 0. 

V x (fjj) = 0. 



(14) 
(15) 



Equations (^J), (|15j) and the boundary conditions are 
satisfied by the Harris equilibrium 



fj = cosh 2 [(xi — 0.5) I a 
J = J e =(0,0, 



B = B< 



P = P = 



tanh(l/2a) cosh 2 [(xi - 0.5) /c 
0i tanh[( !l: 0.5V a ] + --Y 



tanh(l/2a) 



B< 



2 ' 



(16) 
, (17) 

(18) 
(19) 



where i?| and £?f are constants. The resistivity given by 
Eq. (|l^) decreases from the boundary towards the sheet 
center where it takes on a minimum value. This is in 
accordance with the expectation that the plasma is hot- 
ter within the current sheet combined with the decrease 
of the typical Spitzer resistivity with temperature, i.e. 
fj ~ 7 1 - 3 / 2 Unlike other studies |l5|,[l6| where the system 
is infinitely extended in the cross-sheet (xi) direction, 
we do not use the current sheet half width a as the unit 
of length. Instead, we normalize to the finite distance 
L = L\ between the two boundary planes. The magnetic 
field unit, Bq, was chosen in such a way that, in the case 
of .Bl = 0, j-Bf I = 1 on the boundary planes. 
We use the notations 



b = B - B e , j = J — J e , 



(20) 



where v and b are our dynamical variables, for which the 
stress-free boundary conditions are now as follows : 

dv 2 dv 3 db 2 db z 

vi = — - = — - =bi = t— = -— = at xi = 0, 1. 

0x1 0x1 oxi 0x1 

(21) 

We Fourier expand both vector fields into modes ~ 
exp{i(fc2X2 + ^3X3)} in the x 2 and X3 directions. In 
the cross-sheet direction x\ sine and cosine expansions 
are used in correspondence with the imposed stress-free 
boundary conditions (for more details see [B). Dynam- 
ical integrations of the system are performed in Fourier 
space by means of a pseudo-spectral method with 2/3- 
rulc dealiasing. The grid size for the 3D integrations was 
taken to be 32 x 16 x 16 which was found to be sufficient 
for our low Hartmann number studies (see Table 1). Time 
integration was performed using a Runge-Kutta scheme 
with a variable time step. Compared to similar calcula- 
tions with a spatially uniform resistivity, the simulations 
were extremely time-expensive since only very short time 



steps were possible. Additionally one has to keep in mind 
that our spectral resolution is restricted by the evaluation 
of the Jacobian necessary for the linear stability analysis 
of the time-asymptotic states. The used resolution results 
in the inversion of a 2836 x 2836 matrix. 



III. RESULTS 
A. 2D time asymptotics 

We started with a determination of the stability 
boundary for the static sheet pinch equilibrium. The 
Squire theorem allowed a restriction to £3 invariant per- 
turbations (i.e., to perturbations with wave number k% = 
0). Furthermore, due to the invariance of the equilibrium 
in the x 2 direction stability could be tested for each wave 
number k 2 (or the corresponding Fourier mode) sepa- 
rately. First, the system was assumed to be infinitely ex- 
tended in the x 2 direction. In this case the wave number 
k 2 of a perturbation can adopt any real value. Figure |l| 
shows the stability boundary in the Hartmann number- 
wavelength plane for _B| = and a = 0.1 (£?f and a are 
fixed to these values throughout the paper). Since the 
equilibrium profile Bf(xi) is symmetric (i?| = 0), the 
value of i?| has no influence on the stability (see Ref. 
fTo|). The unstable region lies to the right of the bound- 
ary curve. For the spatial resolution used, instability sets 
in at Ha = Ha c = 64.57 and k 2 = k 2c = 2.67. 

In calculations using the full nonlinear equations the 
aspect ratio L 2 (in the 3D case L3 as well) has to be fixed 
to a finite value. We have used L 2 =4 in all nonlinear cal- 
culations. There are some subtleties concerning the onset 
of instability and the application of the Squire theorem 
when L 2 is finite, due to the fact that only a discrete set 
of k 2 values is admitted. With L 2 = 4 instability sets in 
at Ha = [Ha c ]L 2= 4 = 66.20784 and k 2 = [k 2c \L 2 =i = 7r, 
corresponding to a critical wavelength of 2. Unstable 3D 
modes at Hartmann numbers close to [Ha c ]L 2 =i are ex- 
cluded if the aspect ratio L3 is finite (for more details see 
Appendix). Figure ^ shows, for Ha = 70, a comparison of 
the growth rate of the most unstable 2D mode, which has 
wavelength 2 in the x 2 direction, with the growth rates 
of the most unstable 3D mode with the same wavelength 
in the x 2 direction and different wavelengths in the X3 
direction. 

When Ha exceeds the critical value [Ha c ]L 2= 4 the tear- 
ing mode grows due to a bifurcation where a pair of iden- 
tical real eigenvalues becomes positive. A superposition of 
the static equilibrium and the most unstable eigenvector 
was taken as the initial state to follow up the nonlinear 
development of the tearing mode in two spatial dimen- 
sions. 

After several hundred Alfven times ta convergence to 
a stationary state was observed. This state is of course 
linearly stable with respect to two-dimensional pertur- 
bations. It is also clear that the development to the new 
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time-asymptotic states is decelerated the closer to the 
critical value [Ha c ]L 2 =i the Hartmann number is taken. 
This was indicated first by the convergence of the maxi- 
mum eigenvalue to zero. Namely, due to the marginal sta- 
bility with respect to translations in the xi direction, one 
eigenvalue of the time-asymptotic state has to vanish. For 
this real eigenvalue Ao we had, for instance, Ao — — 5T0 -4 
for £ = 70G and #o = 66.21, A ~ -5-1CT 4 for £ = 500 and 
#o = 66.3 (and A ~ -2- 10~ 4 for i = 800 and #0 = 66.3), 
but already A ~ -2 • 10~ 6 for t = 500 and #0 = 67 (and 
A ~ -10~ 8 for i = 800 and Ha = 67). The time develop- 
ment of run 1 became extremely slow. This run close to 
[#a c ]i 2= 4 was performed in order to make it as sure as 
possible that secondary bifurcations close to the primary 
bifurcation point were not overlooked. Since the time- 
asymptotic solutions for Ha = 66.21, 66.3, and 67 are all 
of the same type, the solutions simulated for Ha = 66.3 
and Ha = 67 arc likely to belong to a branch originating 
in the primary bifurcation. In Fig. |^ the time develop- 
ments of the specific kinetic energy E^ in = ^y J v v 2 d V, 
the specific magnetic energy E mag = J v h 2 dV, and 
their sum, the total energy E to t, are plotted for runs 2 
and 3 (run 3 only shown in the inset). Nearly perfect 
steady states are reached for both Hartmann numbers 
at later stages. For Ha = 67 it is practically constant in 
time at t = 800. The amplitude of the eigenvectors is not 
determined by the stability analysis. For Ha = 67 thus 
two energetically different initial conditions were consid- 
ered and were found to relax toward the same asymp- 
totic state, one from energetically below and the other 
from energetically above the asymptotic energy (in the 
first case, not shown in the figure, the energies increase 
as functions of time and then become almost constant). 

In Fig. [| the new asymptotic state is shown for Ha = 
S = M = 67. Field lines of B, stream lines of v, and 
contour lines of the current density component J3 are 
drawn. One observes a magnetic island structure with a 
chain of X and O points, fluid motion in the form of 
convection-like cells or rolls, and a filamentation of the 
original current sheet. For J3 only the most inner part of 
the sheet is shown to highlight the filamentation despite 
of the dominant J| . Two wavelengths in the xi direction 
are seen, corresponding to the fact that L2 = 4 and the 
critical perturbation has wavelength 2. 



B. 3D secondary instability of the 2D 
time-asymptotic states 

Once the 2D time-asymptotic states close to the bi- 
furcation point were calculated with sufficient accuracy, 
their linear stability with respect to 3D perturbations 
could be investigated. Though the stability boundary 
of the quiescent basic state is determined by the Hart- 
mann number, the bifurcating states may depend on S 
and M separately. We have restricted ourselves, how- 
ever, to cases with Ha = S = M . The two-dimensional 



states were extrapolated to three dimensions by contin- 
uing them constantly in the X3 direction, and the stabil- 
ity analysis was performed for the resulting 3D systems. 
Since the equilibria were invariant in the X3 direction, sta- 
bility could be tested for each wave number k% separately. 
Taking into account just one wave number, k$ = =p27r/L3, 
in the x% direction, the aspect ratio or pinch height L3 
was varied. We also added constant magnetic shear com- 
ponents Bf to the saturated 2D states. This could be 
done since such constant field components do not influ- 
ence a 2D solution: The contribution of £>| to the Lorentz 
force JxB in Eq. ([!]) vanishes since both J and the added 
magnetic field are in the £3 direction, and its contribution 
to the term Vx (vxB) in Eq. (g) vanishes as a conse- 
quence of the incompressibility condition V-v = [in the 
incompressible case one has Vx(vxB) = (B-V)v— (v-V)B]. 

The motivation for adding a _B| is that in many appli- 
cations externally generated magnetic fields are present 
in addition to the self-consistently supported ones. In 
the solar atmosphere, for instance, current sheets may 
form when regions of obliquely directed magnetic field are 
brought together and will then in general have a sheet- 
wise field component. In magnetic fusion devices like the 
tokamak toroidal magnetic fields, which correspond to 
sheetwise fields in plane geometry, are externally applied 
to stabilize the confined plasma. 

Results of the stability calculations for Ha = M = S = 
66.3 and Ha = M = S = 67 are shown in Fig. || where the 
maximum real part of the eigenvalue spectrum is plotted 
against the varying parameters L 3 and £?| . In the case of 
S| = the two-dimensional saturated states are always 
unstable, namely to three-dimensional disturbances with 
a sufficiently large wavelength in the £3 direction (see 
upper panel in Fig. |^). 

At the stability threshold always two identical real 
eigenvalues become positive. The multiplicity two results 
from the symmetry of the system with respect to reflec- 
tions in the planes X3 = const., due to which the linearly 
independent modes with wave numbers +£3 and —Aft, re- 
spectively, become simultaneously unstable (the periodic 
boundary conditions, which allow the decomposition into 
Fourier modes, are also needed here) The secondary 
instability to three-dimensional perturbations is always 
suppressed by a sufficiently strong field _B| (see lower 
panel in Fig. |5|). This is in accordance with the general 
expectation that a magnetic field impedes motions with 
gradients in the direction of the field due to the tension 
associated with the lines of force. 

The closer to the critical value [Ha c ]L 2 =4 the Hart- 
mann number is, the larger is the minimum wavelength 
of the unstable perturbations in the third dimension (see 
upper panel in Fig. ^J). Now the Squire theorem does not 
exclude that 3D perturbations to the quiescent basic state 
are unstable immediately above the critical Hartmann 
number, provided their wavelengths 2*k /k% are sufficiently 
large (cf. Appendix). One might suspect, therefore, that 
the unstable 3D perturbations to the 2D time-asymptotic 
tearing-mode state are also unstable perturbations with 
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respect to the basic state at the same Hartmann num- 
ber. This is not the case, however: Consider, for exam- 
ple, the curve for Ha~ 66.3 in Fig. ^ (maximum growth 
rate over wavelength of the perturbation in the X3 di- 
rection). For L3 = 27r/fc3 = 7 one observes 3D instabil- 
ity of the 2D time-asymptotic state. Can the quiescent 
basic state for Ha = 66.3 be unstable to a 3D pertur- 
bation with wavelength 7 in the X3 direction? The wave 
number &2 of such a 3D perturbation can take on the 
values n ■ 2ir/4, n = 1, 2, 3, . . . (since we have chosen the 
fixed aspect ratio £2 = 4). Squire's theorem connects the 
3D perturbation to a 2D perturbation with wave number 
k% = \k\ + (27T/7) 2 ] 1 / 2 which is simultaneously unstable or 
stable at the Hartmann number Ha = (^2/^2)66. 3. With 
fc2 = 27r/4, i.e., with the smallest possible I&2I, one finds 
Ha = 57.6, which is below the critical value Ha c = 64.57 
(see Fig. Q). That is to say, a 3D mode with fc 2 = 27r/4 
(and &3 = 27r/7) cannot be an unstable perturbation to 
the quiescent basic state at _£Za = 66.3. For the next pos- 
sible &2 value, 2 • 27r/4, one has Ha = 63.7, still below 
the critical value Ha c . For all higher ki values the wave- 
length l-Kjki lies clearly below the unstable wavelength 
domain (see Fig. for Ha c <Ha<66.3. 

Figures ^ and |?] show an unstable 3D eigenstate to the 
time-asymptotic 2D state with a = 0.1 and L2 = 4 at 
Ha = 67 (which is shown in Fig. ^). The fields are shown 
in Fig. |6| in the X2-X3 plane to underline qualitatively 
new structures in the third dimension. As in Fig. |J, two 
wavelengths of the perturbation in X2 are shown. The 
2D equilibrium is mixed with the 3D perturbation in the 
ratio 50% equilibrium to 50% perturbation. In the per- 
turbed state velocity and magnetic field have components 
in the X3 direction and all structures, including the cur- 
rent filaments, are modulated in this direction. 

Previous analyses of secondary instabilities of the sheet 
pinch Jl5|,[l6| , as well as analyses of similar instabilities in 
hydrodynamic shear flows jl9| , have indicated that these 
instabilities are ideal, i.e., their growth rates independent 
of dissipation. It appears interesting, therefore, to com- 
pare the growth rates of the secondary instability with 
those of the primary one at the same Hartmann numbers. 
The growth rate of the most unstable 2D tearing mode 
(those with wavelength 2) for iJa = 66.3 is 8.8 • 10~ 4 , the 
corresponding growth rate for Ha = 67 is 7.5 • 10 -3 . A 
comparison with the upper panel in Fig. || shows that 
the secondary instability grows approximately five times 
as fast as the primary one. This agrees with results of 
Dahlburg et al. @ , |l(| . However, since all our calculations 
were restricted to S and R values close to the primary 
bifurcation point, where the growth rates of all primary 
or secondary modes go through zero or are still negative, 
they do not allow yet a characterization of the secondary 
mode as ideal or non-ideal; saturation of the growth rate 
may occur for larger S and R. 



C. 3D time asymptotics 

Finally, full three-dimensional simulations were per- 
formed to follow the unstable modes in their nonlinear 
evolution. The resistivity gradients made the simulations 
again extremely time-expensive. The calculations were 
thus restricted to the case i 3 = L 2 = 4, Bf = _Bf = 0, and 
Ha = 67. The asymptotic 2D state was extrapolated con- 
stantly into the third dimension and was mixed with the 
most unstable 3D eigenstate giving the initial condition. 
The first phase of the full three-dimensional simulation 
was done with a lower spectral resolution, namely 16 3 , up 
to the time to ~ 2900. In Fig. ||, the temporal behavior of 
the specific energies after this initial growth phase for the 
next 900 time units is shown which was calculated with 
the highest resolution as given in Table |. The energies 
still oscillate slightly, but with a decreasing amplitude, in- 
dicating convergence to a three-dimensional steady state 
of the sheet pinch configuration. 

Furthermore, the solution is characterized by a clear 
and apparently time-independent spatial structure. In 
Fig. |l^ corresponding level surfaces of |v| are shown at 
t = 550. We found the same shape of the level surface of 
|v| at t = 900. Additionally, the same level surfaces are 
shown for the 2D time-asymptotic state with otherwise 
the same parameters in Fig. ^|. The comparison indicates 
that there is some relation between the two solutions — 
the isosurfaces of |v| in the 3D case are obtained from 
those in the 2D case by a modulation in the X3 direction. 
This suggests, but does not prove, that the unstable 3D 
perturbations to the time-asymptotic 2D state do not 
drive the system to a completely different solution ex- 
isting somewhere in phase space, but that 2D and 3D 
solutions originate simultaneously in the primary bifur- 
cation of the quiescent basic state. 



IV. SUMMARY AND OUTLOOK 

We have numerically studied the primary and sec- 
ondary bifurcations of an electrically driven plane sheet 
pinch with stress- free boundaries. The profile of the elec- 
trical conductivity across the sheet was chosen such as to 
concentrate the electric current largely about the mid- 
plane of the sheet and thus to allow instability of the 
quiescent basic state at some critical Hartmann number. 
Our results can be summarized as follows. 

(1) The most unstable perturbations to the basic state 
are two-dimensional tearing modes. If the whole prob- 
lem is restricted to two spatial dimensions, also the bi- 
furcating new time-asymptotic state is of the tearing- 
mode type, namely, a stationary solution characterized 
by a magnetic island structure with a chain of X and O 
points, fluid motion in the form of convection-like rolls, 
and a filamentation of the original current sheet. We have 
calculated this state with precision for the aspect ratio 
L2 = 4 and Hartmann numbers close to the critical one. 
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In contrast to the stability boundary the bifurcating new 
solutions may depend on the two Reynolds-like numbers 
of the problem separately. We have restricted ourselves 
to M =S= Ha. 

(2) The bifurcating steady (time-asymptotic) two- 
dimensional state was tested for stability with respect to 
three-dimensional perturbations. It proved to be unstable 
to 3D perturbations with a sufficiently large wavelength 
in the third direction. At the stability threshold always 
two identical real eigenvalues become positive (i.e., there 
are two purely growing unstable eigenmodes). We also 
added constant external magnetic field components along 
the invariant direction to the 2D tearing-mode equilib- 
rium. If these components are sufficiently strong, they 
suppress the secondary instability with respect to three- 
dimensional perturbations, which is in accordance with 
the general expectation that a magnetic field impedes 
motions with gradients in the direction of the field and 
has been noted before jl6) . 

(3) Full three-dimensional simulations were performed 
to follow the unstable 3D modes in their nonlinear evolu- 
tion. The solution seems to converge to a 3D steady state. 
Although velocity and magnetic field have now compo- 
nents in the invariant direction of the 2D state and all 
structures are modulated in this direction, there is still 
some resemblance to the 2D tearing mode state. This 
suggests, but does not prove, that the unstable 3D per- 
turbations to the 2D state do not drive the system to a 
completely different region in phase space. 2D and 3D 
solutions might originate simultaneously in the primary 
bifurcation of the basic equilibrium. 

Since our calculations were made very close to the 
primary bifurcation point, we suppose that the steady 
2D tearing-mode state is unstable from the beginning 
and that there is a direct transition of the system from 
the quiescent basic state to three-dimensional attractors. 
This is true if the 2D state is not stabilized by an external 
magnetic field along its invariant direction. Furthermore, 
a sufficiently small L$ ensures stability of the 2D state 
in a certain Hartmann number interval above the critical 
value (since the unstable 3D perturbations, whose wave- 
length must exceed some threshold value, are then not 
admitted). Finally, the aspect ratio L2 and the magnetic 
Prandtl number Pr m = u/rio = S/M might influence 
the bifurcation scenario, possibly in such a way that in 
some parameter ranges the 2D tearing-mode solution bi- 
furcates stably from the basic state. 

Secondary instabilities that succeed primary two- 
dimensional ones and that lead to three-dimensionality 
have been considered as an important step in the transi- 
tion from laminar to turbulent states in linearly unstable 
nonconducting shear flows []l9|-^f||. For the case of lin- 
early stable shear flows (e.g. the plane Couette flow), 
it was suggested that nonlinear stationary and linearly 
unstable three-dimensional states, which develop already 
below the onset threshold of turbulence |22|,^3| , can form 
a chaotic repellor in phase space |24[ ]. Such a repellor 
can cause the transient turbulent states above the on- 



set threshold. There is some analogy of the magnetohy- 
drodynamic pinch to shear flows, and Dahlburg et al. 
]f5| , [T6[ have presented numerical evidence that the sec- 
ondary instability of two-dimensional quasi-cquilibria of 
the tearing-mode type can lead to turbulence in a plane 
sheet pinch. In their calculations the pinch was not driven 
by an external electric field (nor mechanically driven) and 
the electrical conductivity was assumed to be spatially 
uniform. In such a case the pinch always decays resis- 
tively, that is, velocity and magnetic field tend to zero 
as t — > 00. By our choice of the resistivity profile and 
the applied electric field we could calculate exact time- 
asymptotic, in particular steady states and could corrob- 
orate the result of Dahlburg et al. that saturated two- 
dimensional tearing-mode states are unstable to three- 
dimensional perturbations. We did not observe a tran- 
sition to a turbulence-like state yet. Irregular behavior 
may be expected to arise through subsequent bifurca- 
tions when the Reynolds-like numbers are further raised. 
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APPENDIX: INSTABILITY OF THE QUIESCENT 
BASIC STATE AND SQUIRE'S THEOREM IN 
THE CASE OF A FINITE ASPECT RATIO L 2 

Squire's theorem states that for increasing Hartmann 
number two-dimensional perturbations to the quiescent 
basic state become unstable first. Specifically: For each 
three-dimensional eigenmode with wave numbers &2, ky, 
and growth rate A at Hartmann number Ha, there ex- 
ists a two-dimensional (i.e., X3 invariant) eigenmode with 
wave number &2 = (&f + ^I) 1 ^ 2 an( ^ growth rate A = 
(k~2/k2)\ at Hartmann number Ha—(k,2/k2)Ha B. 

1. Case L2 — 00 

If the stability problem is considered on the infinite X2- 
X3 plane, i.e., with all wave numbers and k^ allowed, 
for increasing Ha one or several two-dimensional modes 
with a critical wave number k2 c (and the correspond- 
ing critical wavelength L2 C = 27r//c2 c ) become unsta- 
ble at a critical Hartmann number Ha c where all three- 
dimensional modes are still stable. Above Ha c the critical 
value &2c broadens to an unstable k2 interval. The latter 
means, however, that three-dimensional modes could be 
unstable immediately above Ha c . Namely, consider a 3D 
eigenmode with wave numbers k2, k% and growth rate 
A at some Hartmann number Ha — Ha c + e, e > 0. If 
/c2 is chosen from the interior of the unstable k2 interval 



G 



at Ha, then fe = (fcf + k 2 ) 1 / 2 lies within the unstable 
fe interval at the Hartmann number Ha — (fe/fe)iJa, 
where Ha c < Ha < Ha c + e, if only |fe| is chosen suffi- 
ciently small. This does not mean yet that the 2D mode 
to which the 3D mode is connected is unstable, since 
there are in general also stable 2D eigenmodes with the 
same wave number fe . But if the associated 2D mode is 
unstable, i.e., if the real part of A= (fe/fe)A is positive, 
this implies that also 5R(A) > 0. The possibility of un- 
stable three-dimensional eigenmodes close to the critical 
Hartmann number is excluded by the Squire theorem, 
however, if L3 is finite, i.e., if there is a positive lower 
bound (however small) to the modulus of the wave num- 
ber fe. In that case there is a finite Hartmann number 
interval above Ha c where all unstable eigensolutions are 
purely two-dimensional. 



2. Case L2 finite 

Fixing Li to a finite value complicates the problem, 
since only a set of discrete values is admitted for fe. If 
not just L2 = n-27r/fe c , with n denoting a positive integer 
number, that is, if fe c is not just an admissible fe, insta- 
bility to 2D modes will set in at some Hartmann number 
[Ha c ]L 2 above Ha c and for a wave number [fe c ]_L 2 differ- 
ent from fe c . 



a. Subcase L2 < L2 C = 2ir/k2 c 

In the case L2 < -^2 C I0r a U fe holds fe > fe c and 
consequently the smallest admissible fe becomes unsta- 
ble first, i.e., [fe c ]_L 2 — 2tt/L2 > fe c . It is easily seen that 
as in the case of L2 = 00 (i) directly at the onset of insta- 
bility only 2D modes can be unstable (since modes with 
fe = cannot be unstable JTo| ] and the Squire theorem 
thus would connect any unstable 3D mode to a 2D mode 
with wave number fe > [fe c ].L 2 outside the unstable fe 
interval at the Hartmann number [Ha c ]L 2 ), (ii) immedi- 
ately above [Ha c ]L 2 also unstable 3D modes are possible 
(or at least not forbidden by Squire's theorem), (iii) a 
finite aspect ratio L3 (however large) ensures that in a 
finite Hartmann number interval close to the onset of 
instability only purely two-dimensional eigenmodes are 

unstable (see also Fig. where the pinch is stable with 

< 



respect to 3D modes for L3 < 6). 



b. Subcase L2 > L2 



27v/k 2 



More involved is the situation for L2 > L2 C . Then 
it cannot be excluded generally that 3D modes become 
unstable first, and in principal each individual situation 
has to be tested separately. One can distinguish between 
the cases [fe c ]_L 2 > fe c and [fe c ]_L 2 < fe C j of which 



the first one is simpler. In both cases special compli- 
cations arise from the fact that 3D modes with wave 
numbers fe smaller than [fe c ]z, 2 , that is to say, with 
k 2 = n ■ 2tt/L 2 < [fe c ]L 2 = n o ■ 2tt/L 2 (n, n denoting 
integer numbers) can come into play. 

In the case of [fe c ]L 2 > fe c these 3D modes (with 
k2 smaller than [fe c ]z, 2 ) are the only 3D modes that 
could become unstable at a Hartmann number less than 
[Ha c }L 2 (where the first 2D mode becomes unstable); 
if th ey rem ained stable, the situation is similar to sub- 



case 



A2; 



The 3D modes with fe < fe c must remain 
stable close to the onset of 2D instability, however, if 
[Ha c }L 2 does not exceed Ha c too much (and [fc2 c ]i 2 does 
not differ too much from £)2 C ), such that (i) \k%\ has to be 
larger than some positive threshold value in order that 
k2 = [k\ + ki) 1 / 2 (with k2=n ■ 2n/ L2, n < no) can come 
into the unstable fe interval close to the onset of insta- 
bility (since there is a finite gap between the unstable k% 
interval and the largest admissible fe that is smaller than 
k2 c ) and (ii) as a consequence of this Ha — (fe/A^i/a 
must be smaller than Ha c . If this is the case and, fur- 
thermore, [fc2 c ]L 2 > fee: th e situation is the same as for 
L 2 < 2ir/k 2c . 

The numerical example of this paper belongs to the 
category just discussed: L2 finite, L2 > 2ir/k 2c , [feJi-2 > 
fe c , and close to the onset of instability no unstable 3D 
modes with fe < [fe c ]z, 2 . We found fe c = 2.67, cor- 
responding to a critical wavelength of L 2c = 2.35, and 
Ha c = 64.57. The critical values for the fixed aspect 
ratio L2 = 4 are [fe c ]z, 2 =4 = i", corresponding to a 
critical wavelength of [i2 c ]i 2 =4 = 2, and [Ha c ]L 2 =4 = 
66.20784 (see also Fig. [I]). Loosely speaking, an unsta- 
ble 3D mode has to fit now between Ha c and [-ffa c ]i 2= 4 
with its critical Hartmann number Ha. It can only have 
the wave number fe = 27r/4 = 7r/2, since otherwise 
fe = (fc| + fcf) 1 / 2 > [fe c ]i 2 (he., the associated 2D mode 
would be stable). This implies, in order to have instability 



:[(7T/2) 



and consequently 



fe 2 ] 1/2 >fe 



= 2.67 , 



(Al) 



Ha=^[Ha c ] Lo=i < J— [Ha c i 2=4 ~38.9 , (A2) 

fe " fee 



which lies below Ha c . 3D modes with fe = tt/2 can- 
not be unstable even for Hartmann numbers signifi- 
cantly above [-Ho c ]l 2= 4. 3D modes with fe = [fe c ]L 2 =4j 
on the other hand, can be unstable immediately above 
[Ha c ]L 2= 4 and are stabilized by an upper bound to the 
aspect ratio L3 , as discussed in the preceding subsections 
of this Appendix. We found the condition L3 < 1000 to 
be sufficient to stabilize all 3D modes at Ha = 66.208 
(> [ffa c ]z, 2=4 = 66.20784). 
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TABLE I. Parameters of the 2D and 3D sheet pinch simu- 
lations. 



Run 


Ha 


L 2 




(Vi,A 2 ,A 3 ) 


1 


66.21 


4 




(32,16,-) 


2 


66.30 


4 




(32,16,-) 


3 


67.00 


4 




(32,16-) 


4 


67.00 


4 


4 


(32,16,16) 
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FIG. 1. Stability boundary for the quiescent basic state 
with a — 0.1 and _B| = in the Ha-L-2 plane. L2 — 2n/k2 is the 
wavelength of the perturbation in the X2 direction (stability 
was tested for each wave number &2 separately). 
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FIG. 3. Temporal behavior of the specific kinetic (dotted 
line), magnetic (dashed line), and total (solid line) energies 
for the nonlinear two-dimensional development of the tear- 
ing mode for Ha = S = M — 66.3 , a = 0.1, and L 2 = 4. The 
inset shows additionally the corresponding development for 
Ha = S=M = &7 . 
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FIG. 2. Maximum real parts of the eigenvalue spectrum 
vs. Li resulting from the stability analysis of the quiescent 
basic state for a — 0.1 and L2 = 2 at Ha = R = S = 70. 
L3 = 2n/k3 is the wavelength of the perturbation in the X3 
direction (stability was tested for each wave number fc 3 sepa- 
rately). The dotted line marks the maximum growth rate of 
the most unstable 2D mode (ks — 0) for the same parameters. 



FIG. 4. Magnetic field lines (left), velocity stream lines 
(middle), and contour lines of the current density component 
J 3 (right) for Ha = 67, a = 0.1, and L 2 = 4 at t = 800. Solid 
(dashed) velocity stream lines correspond to clockwise (coun- 
terclockwise) motion. Only the inner part of the current sheet 
is shown for the current density contour plot to highlight the 
filamentation of J3. 
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FIG. 6. Most unstable 3D eigenstate with wave number 
h A = 2-7T-/4 to the time-asymptotic 2D state (see also Fig. ^ at 
Ha — 67. Vector plots of the velocity field (left) and contour 
lines of the magnetic field component _B 3 (right) in planes 
X\ — 0.6 (upper row) and x\ = 0.4 (lower row) are shown. 
Dotted contour lines indicate negative B3, solid lines positive 
ones. The time-asymptotic 2D state is mixed with the 3D per- 
turbation in the ratio 50% equilibrium to 50% perturbation. 



FIG. 7. (JPG-File) Isosurfaces |J|=5.45 for the same data 
set as in Fig. ^[ The maximum value of |J| is 9.77. 



FIG. 5. Maximum real parts of the eigenvalue spectrum 
versus L 3 (upper panel) and B| (lower panel) for a — 0.1 and 
1/2 = 4 at different Hartmann numbers Ha (always 5 = M). 
1/3 = 2n/ks is the wavelength of the perturbation in the X3 
direction (stability was tested for each wave number £3 sepa- 
rately) . 
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FIG. 8. Temporal behavior of the specific kinetic, mag- 
netic, and total energies for the nonlinear development of the 
three-dimensional system for L2 = L3 = 4, Bf = _B| = 0, and 
Ha = 67. The time to was 2900. 



FIG. 9. (JPG-File) Isosurfaces |v| =0.03 and |v| =0.016 for 
the time-asymptotic 2D state. The values of the parameters 
are L2 = 4, £>| = £>| = 0, and Ha = 67. The maximum and 
minimum values of |v| are 0.0384 and 0.0017, respectively. 

FIG. 10. (JPG-File) Isosurfaces |v|=0.03 and |v|= 0.016 
for the selfconsistent 3D state at t = 550. The values of the 
parameters are L2 = L3 = 4, S| = S| = 0, and Ha = 67. 
The maximum and minimum values of |v| are 0.0311 and 0.0, 
respectively, in the 3D case. 
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